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The nature of QCD phase transition is studied with massless up and down 
quarks and a light strange quark, using the Wilson formalism for quarks on 
a lattice with the temporal direction extension N t = 4. We find that the 
phase transition is first order in the cases of both about 150 MeV and 400 
MeV for the strange quark mass. These results together with those for three 
degenerate quarks suggest that QCD phase transition in nature is first order. 



One of major goals of numerical studies of lattice QCD is to determine 
the nature of the phase transition from the high temperature quark-gluon- 
plasma phase to the low temperature hadron phase, which is supposed to 
occur at the early stage of the Universe and possibly at heavy ion collisions. 
It is, in particular, crucial to know the order of the transition to understand 
the evolution of the Universe. Because the critical temperature is of the 
same order of magnitude as the strange quark mass, we have to take into 
account the effect of the strange quark as well as those of almost massless 
up and down quarks. In this article we study the nature of QCD phase 
transition with these three quarks in the Wilson formalism of fermions on the 
lattice. The investigation using Wilson quarks is highly important, because 
the Wilson formalism of fermions on the lattice is the only known formalism 
which possesses a local action for any number of flavors with arbitrary quark 
masses. Preliminary reports are given in [|TJ. 

The results of numerical simulations with Wilson quarks as well as with 
staggered quarks |§ imply that the chiral phase transition is second order for 
two massless quarks, while it is first order for three. This is also consistent 
with theoretical predictions based on universality |4j]. In a numerical study 
we are able to vary the mass of the strange quark. When the mass of the 
strange quark is reduced from infinity to zero, the phase transition must 
change from second order to first order at some quark mass m*. This is 
a tricritical point. The crucial point is whether the physical strange quark 
mass is larger or smaller than m*. 

The action we use in this article is composed of the standard one-plaquette 
gauge action with the gauge coupling constant g {(3 = 6/g 2 ) and the Wilson 
quark actions for u, d and s quarks with the hopping parameters K u , K d and 
K s , respectively. We set K u = K d . 

We perform simulations on lattices of the temporal extension N t = 4 with 
the spatial sizes 8 2 x 10 and 12 3 . Although Nt = 4 is not large enough to 
obtain the result for the continuum limit, this work is a first step toward the 
understanding of the nature of the QCD phase transition with Wilson quarks. 
We generate gauge configurations by the hybrid R algorithm with a molecular 
dynamics step At = 0.01 and a trajectory of unit time. The inversion of the 
quark matrix is made by a minimal residual method or a conjugate gradient 
method. When the hadron spectrum is calculated, the lattice is duplicated 
in the direction of lattice size 10 or 12. We use an anti-periodic boundary 
condition for quarks in the t direction and periodic boundary conditions 
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otherwise. The statistics is in general total r = several hundreds, and the 
plaquette and the Polyakov loop are measured every simulation time unit 
and hadron spectrum is calculated every 5r = 5 (or 10 sometimes). When 
the value of (3 is small, the fluctuation of physical quantities are small and 
therefore we think the statistics is sufficient for the purpose of this article. 

Let us first discuss the results for the degenerate N F = 3 case: K u = 
Kd = K s = K. In order to find the transition/crossover points we perform 
simulations at /3=4.0, 4.5, 4.7, 5.0 and 5.5. The simulation time history of the 
plaquette at (3 = 4.7 on a 4 x 12 3 lattice is plotted in Fig. [l|(a). The low and 
high temperature phases coexist over 1,000 trajectories at K = 0.1795 and 
in accord with this we find two-state signals also in other observables such as 
the plaquette and m\ (Fig. From them we conclude that the transition 
at K = 0.1795(5) and f3 = 4.7 is first order. The transition/crossover points 
thus identified are summarized in Table 1. 

In a previous paper we reported a clear two-state signal in the chiral 
limit K c = 0.235 at f3 = 3.0 for N F = 3. An interesting question is to 
determine the critical value of the quark mass m c J' lt up to which the first 
order phase transition persists. We observed clear two state signals at (3 = 
4.0, 4.5 and 4.7, while for (3 = 5.0 and 5.5 no such signals have been seen (see 
Fig. |l|(b) for the time history of the plaquette at (3 — 5.0). At the transition 
point (in the confining phase) of (3 = 4.7 the value of m q a is 0.175(2) and 
(m^/uip) = 0.873(6). The quark mass is defined through the axial Ward 
identity || for Wilson quarks. The results of the hadron spectrum in the 
range of (3 — 3.0 - 4.7 for N F = 2 and 3 (Fig. |3|) indicate that the inverse 
lattice spacing a~ l estimated from the p meson mass is almost independent 
on (3 in this range and a -1 ~ 0.8 GeV. Therefore we obtain a bound on the 
critical quark mass m c ^ lt > 140 MeV, or equivalently (m n /m p ) crit > 0.873(6). 

We note that these values are much larger than those for the staggered 
fermion case for which m° rit ~ 10 - 40 MeV (m^ rit a = 0.025 - 0.075 g [7j . We 
have used a^ 1 ~ 0.5 GeV at (3 = 5.2 for V F =2§) and (m^/m p ) crit ~ 0.42 - 
0.58 (the result for N F = 4 at [3 = 5.20, because the data for N F = 3 are 
not available.) 

Now let us discuss a more realistic case of massless u and d quarks and 
a light s quark [Np = 2 + 1). Our strategy to investigate this problem is 
similar to that we utilized in our previous work for the investigation of the 
chiral transition in the degenerate quark mass cases, which we call "on-i^ c " 
simulation method. We set the value of the masses of u and d quarks m u d 
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to zero (K u = Kd = K c ) and fix the s quark mass m s to some value, and 
make simulations reducing the value of (3 starting from a deconfining phase. 
When u and d quarks are massless, the number of iteration 7V inv needed for 
the quark matrix inversion (for u and d quarks) is enormously large in the 
confining phase, while it is of order of several hundreds in the deconfining 
phase. This is due to the fact that there are zero modes around K c in the 
low temperature phase, while none exists in the high temperature phase [|10|1 . 
Therefore the iVi nv provides an extremely good indicator to discriminate the 
two phase when u and d quarks are on the K c line. Combining this with 
measurements of physical observables, we identify the transition point and 
determine the order of the transition. 

The values which we take for K c are given in Table 2. They are the 
vanishing point of extrapolated m 2 for Np = 2 and interpolated ones. We 
have used those for Np = 2, because we have the data most in this case, 
and the difference between that for Np = 3 and for Np = 2 is of the same 
order of magnitude as the difference due to the definition of K c , either the 
vanishing point of m q or m 2 0. 

We study two cases of m s ~ 150 MeV and 400 MeV. From the value of 
a -1 ~ 0.8 GeV and an empirical rule m q a ~ (2/3) (1/ i^T — 1/K C ) satisfied 
in the (3 region we have studied (Fig. ||), we get the values for K s shown 
in Table 2. It should be noted that the physical s quark mass determined 
from = 1020 MeV turns out to be m s ~ 150 MeV in this (3 range for our 
definition of the quark mass. 

The simulation time history of N inv on the 8 2 x 10 spatial lattice is plotted 
in Fig. f|(a) for the case of 150 MeV. When (3 > 3.6, N inv is of order of several 
hundreds, while when (3 < 3.4, N\ nv shows a rapid increase with r. At (3 = 3.5 
we see a clear two-state signal depending on the initial condition: For a hot 
start, iVi nv is quite stable around ~ 900 and m 2 is large (~ 1.0 in lattice 
units). On the other hand, for a mix start, A^ inv shows a rapid increase with 
r and exceeds 2,500 in r ~ 10, and in accord with this, the plaquette and 
m 2 decreases with r as shown in Fig. |](b) for the plaquette. For the case 
of 400 MeV a similar clear two-state signal is observed at (3 = 3.9 both on 
the 8 2 x 10 and 12 3 spatial lattices (Fig. |). Thus the method of "on-K c " 
simulation is very powerful to identify the point of the transition, although 
in the low temperature phase we are only able to obtain upper (or lower) 
bounds of physical observables. The values of m 2 versus j3 are plotted in 
Fig. ^ together with those in the case of degenerate Np = 3 on the K c line 
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0]. At /3 = 3.5 for the case of 150 MeV and at j3 = 3.9 for the case of 400 
MeV, we have two values for m\ depending on the initial configuration. The 
larger ones of order 1.0 are for hot starts, while the smaller ones are upper 
bounds for mix starts. 

It should be noted that the critical (3 value of 3.9 for the case of m s ~ 400 
MeV is close to that of the chiral transition for Np = 2 at (3 ~ 3.9 - 4.0 
0. When m s is larger enough than the inverse lattice spacing, the transition 
point should become almost identical with that of Np = 2 and the transition 
should change from first order to continuous. Therefore the closeness of the 
two critical points is quite natural, because the value of 400 MeV is not far 
from that of the inverse lattice spacing of about 800 MeV. 

The above results imply that 

m* > 400 MeV. 

Our results for the presence or absence of the transition in the cases of de- 
generate Np = 3 and Np = 2 + 1 are summarized in Fig. |7|. Clearly the 
point which corresponds to the physical values of the u,d and s quarks exists 
in the range of the first order transition. Thus our results suggest that the 
transition for the physical quark masses is of first order. 

The result by the Columbia group [0] for staggered quarks shows that no 
transition occurs at m u a = m^a = 0.025, m s a = 0.1 (m u = m<j ~ 12 MeV, 
m s ~ 50 MeV using a -1 ~ 0.5 GeV), which suggests that a first-order phase 
transition does not occur for the physical u,d and s quarks. This result is in 
conflict with ours, albeit with slightly heavier masses for u and d quarks. One 
possibility for the discrepancy is that both/either are far from the continuum 
limit. Consistency between Wilson and staggered quarks is an issue which 
should be investigated in future. 

The simulations are performed with HITAC S820/80 at KEK and with 
VPP500/30 at the University of Tsukuba. We would like to thank members 
of KEK for their hospitality and strong support. This project is in part 
supported by the Grant-in- Aid of Ministry of Education, Science and Culture 
(No.06NP0401). 
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p 


K t 


3.0 


> 0.230 


4.0 


0.200 - 0.205 


4.5 


0.186 - 0.189 


4.5* 


0.186 - 0.189 


4.7* 


0.179 - 0.180 


5.0 


0.166 - 0.167 


5.0* 


0.166 - 0.1665 


5.5 


0.1275 - 0.130 



Table 1: Finite temperature transition/crossover K t for Np = 3 obtained on 
an 8 2 x 10 x 4 lattice (data with * obtained on a 12 3 x 4 lattice). At (3 < 4.7 
clear two-state signals are observed on both lattices. 



m s « 150MeV 
P K ud K s 


m s ~ 400MeV 

P K ud K s 


3.2 0.2329 0.2043 

3.4 0.2306 0.2026 

3.5 0.2295 0.2017 

3.6 0.2281 0.2006 
4.0 0.2226 0.1964 


3.7 0.2267 0.1692 

3.8 0.2254 0.1684 

3.9 0.2240 0.1677 
4.0 0.2226 0.1669 
4.3 0.2180 0.1643 



Table 2: Hopping parameters for N F = 2 + 1 simulations performed on 
8 2 x 10 x 4 and 12 3 x 4 lattices. K ud for u and d quarks is set to be equal to 
K c and K s for s quark is chosen so that m s w 150 MeV and 400 MeV in the 
left and right columns, respectively. 
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Figure 1: Time history of the plaquette for Np = 3 at (a) (3 — 4.7 and (b) 
5.0 on a 12 3 x 4 lattice. 
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Figure 2: (m^a) 2 and the plaquette for N F — 3 at f3 — 4.7. 
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Figure 3: m p a and 2m q a in the low temperature phase. Open symbols are 
for N F = 2, = 3.0, 3.5, 4.0, 4.3, and 4.5 on an 8 2 x 10 x 4 lattice. Filled 
symbols are for N F = 3, f3 = 4.0, 4.5 and 4.7 on 8 2 x 10 x 4 and 12 3 x 4 
lattices. K c for Np = 2 is used. Horizontal errors are from those of K c with 
taking into account the difference due to definitions, either the vanishing 
point of m 2 or m q . 
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Figure 4: Time history of (a) A^ inv and (b) the plaquette for m s ~ 150 MeV 
on an 8 2 x 10 x 4 lattice. 
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Figure 5: Time history of the plaquette for m s ~ 400 MeV on a 12 3 x 4 
lattice. 
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Figure 6: (m n a) 2 versus (3 for m s ~ 0, 150 and 400 MeV with m u d ~ 0. 
Filled and open symbols are for 8 2 x 10 x 4 and 12 3 x 4 lattices, respectively. 
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N p = 2 pure gauge 




Figure 7: Order of the finite temperature QCD transition in the (m u d, m s ) 
plane. First order signals are observed at the points marked with filled circle, 
while no clear two state signals are found at the points with open circle. 
The values of quark mass in physical units are computed using a -1 ~ 0.8 
GeV. The real world corresponds to the point marked with star. The second 
order transition line is suggested to deviate from the vertical axis as m U( i oc 
(m* — m s ) 5 / 2 below to*. 
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